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I.  INTRODUCTION 


The  primary  motivation  for  this  work  is  to  explore  turbulence  models  for 
high  speed,  three-dimensional  flows  such  as  those  which  occur  in  the  vicinity 
of  artillery  projectiles.  The  basic  geometry  is  axisymmetric  but  at  angle  of 
attack  or  yaw.  At  high  spin  rates,  the  three-dimensional  effects  are 
significant.  The  range  of  angles  of  attack  is  usually  small  for  a  stable 
projectile  but  practical  geometries  often  involve  surface  irregularities  such 
as  fuses  and  rotating  bands.  These  surface  perturbations  sometimes  lead  to 
local  regions  of  separation  so  that  it  is  desirable  that  the  model  of 
turbulence  be  extendable  to  such  conditions.  Finally  the  projectile  base  flow 
has  a  very  significant  effect  on  the  projectile  drag  and  therefore  a  generally 
useful  model  should  be  applicable  in  the  wake  region. 

The  turbulence  model  must  be  applicable  to  compressible  flow  because  most 
artillery  launching  speeds  are  supersonic.  The  lower  end  of  the  speed  range, 
however,  must  also  encompass  transonic  Mach  numbers  because  this  is  often  a 
critical  flight  condition  which  must  be  studied  in  detail. 

Numerical  techniques  for  investigating  projectile  flow  fields  have  been 
developing  rapidly  over  the  past  15  years;  starting  with  Euler  inviscid  flow 
solutions  coupled  with  boundary  layer  techniques  and  currently  concentrating 
on  various  forms  of  the  thin  shear  layer  Navier- Stokes  codes  which  solve  the 
inviscid  and  viscous  field  at  the  same  time.  The  methods  of  primary  interest 
for  the  turbulence  modeling  considered  here  are:  the  three-dimensional 
boundary  layer  technique  developed  by  Dwyer  and  Sanders,1  the  parabolized 
Navier-Stokes  technique  of  Schiff  and  Steger2  with  the  numerical  algorithm 
based  on  the  work  of  Beam  and  Warming3  as  well  as  the  unsteady  transonic 
Navier-Stokes  technique  of  Pulliam,  Steger  and  Nietubicz.4  5 


2.  hcyer,  H.A.  ,  and  Sanders,  B.E.,  "Magnus  Faroes  on  spinning  supersonic 
Cones.  Part  1:  the  Boundary  Layer,"  BRL  Contractor  Report ,  ARBRL-C2-248 , 
July  1975 . 

2.  Schiff,  L.B.,  and  Steger ,  J.L.,  " Numerical  Simulation  of  Steady  Supersonic 
Viscous  Flop,"  AIAA  Paper  79-0130 ,  January  1979 . 

3.  Beam,  R.M. ,  and  Warming,  R.F.,  ", An  Implicit  Factored  Scheme  for  Hyperbolic 
Systems  in  Conservation-Lap  Form."  J.  of  Computational  Phusics.  Vol.  22, 
1976. 

4.  Pulliam,  T.H.,  and  Steger,  J.L.,  "On  Implicit  Finite-Difference  Simulations 
of  Three-Dimensional  Flop."  AIAA  Journal .  Vol.  18,  2,  pp.  159-167 , 

February  1980. 

6.  Hietubicz,  C.J.,  Pulliam,  T.H.,  and  Steger,  J.L.,  "numerical  Solution  of  the 
Azimuthal-Invariant  Thin-Layer  Navier-Stokes  Equations ,"  AIAA  Paper  79- 
0010 ,  January  1979. 
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Initially  these  numerical  methods  have  employed  algebraic  turbulence 
models  based  on  the  mixing  length  theory  of  Prandtl.6  These  have  been  modi¬ 
fied  for  compressibility  and  extended  to  include  effects  of  heat  transfer, 
pressure  gradient  and  many  other  effects.  One  of  the  main  criticisms  of  these 
methods  is  that  they  are  local  models  which  do  not  account  for  the  convection 
and  diffusion  of  turbulence.  In  order  to  include  more  of  the  physics  of  tur¬ 
bulence  in  the  mathematical  models,  the  1-  and  2-differential  equation  tech¬ 
niques  were  developed.  A  number  of  surveys7  8  9  of  the  development  and  status 
of  all  these  methods  are  available  and  therefore  no  discussion  is  required 
here.  The  K-E  model  of  Jones  and  Launder10  was  selected  as  the  non-local 
model  to  be  tested  in  the  above  numerical  codes  as  a  representative  differ¬ 
ential  equation  technique. 

This  report  consists  of  three  main  sections:  the  first  contains  a  review 
of  the  three-dimensional,  turbulent  equations  of  motion  aimed  at  clarifying 
the  assumptions  inherent  in  their  subsequent  use,  the  application  to  the 
three-dimensional  boundary  layer  equations  and  finally  their  application  to 
the  steady,  parabolized  Navier-Stokes  equations.  The  inclusion  of  the  K-E 
equation  to  the  boundary  layer  is  based  on  the  work  of  van  Gulick11  12  and  the 
application  to  the  PNS  equations  has  been  reported  by  Kim. 13  Concurrently, 


e.  Schlichting ,  H.  ,  "Boundary  Layer  Theory,"  7  th  Edition ,  McGraw-Hi 1 1  Book 
Co.,  19?9. 

7.  Launder,  E.E.,  and  Spalding,  D.B.,  "Mathematical  Models  of  Turouler.ee," 
Academic  Press,  Inc.,  1972. 

8.  Pxbesin,  K.  W. ,  "numerical  Turbulence  Modeling ,"  A  GAUD  Lecture  Series  Vc.ti' 
on  Computational  Fluid  Dynamics ,  pp.  3-1  to  3-37 ,  April  1977, 

9.  Reynolds,  V.C.,  "Computation  cf  Turbulent  Flows,"  Ann.  Rev.  Fluid  Me ch. , 

Vol.  8,  pp.  183-208 ,  1976. 

10.  Jones,  W.F.,  and  Launder,  B.E. ,  "The  Calculation  of  Lou  Reynolds  Lur.ter 
Phenomena  with  a  Two  Equation  Model  of  Turbulence ,"  Int.  J.  Heat  and  Mass 
Transfer,  Vol.  16,  1973. 

11.  Van  Gulick ,  P. ,  "Application  of  the  K-E  Turbulence  Model  to  a  Boundary 
Layer  Solution  for  Flow  About  a  Spinning  Jawed  Projectile  at  Mach  3 ,” 
Masters  Thesis,  Mechanical  and  Aerospace  Engineering  Department,  Univer¬ 
sity  of  Delaware,  June  1983. 

12.  Van  Gulick,  P. ,  and  Dahberg ,  J.E.,  "Application  of  the  K-E  Turbulence  Model 
to  a  Boundary  Layer  Solution  for  Flow  About  a  Spinning  Jawed  Projectile 

at  Mach  3,"  Proceedings  12th  Southeastern  Conference  on  Theoretical  and 
Applied  Mechanics ,  May  1984. 

13.  Kirri*  vT’»  " Computation  of  Three-Dimensional  Turbulent  Flow  with  Parabolized 
Navier-Stokes  Equations  and  K-E  Turbulence  Model,"  Masters  Thesis, 
Mechanical  and  Aerospace  Engineering  Department,  University  of  Delaware, 
January  1984. 
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the  K-E  technique  has  been  tested  in  the  unsteady,  transonic  Navier-Stokes 
problem  by  Sahu*4  but  discussion  of  his  work  is  not  considered  here. 


II.  TURBULENT  COMPRESSIBLE  EQUATIONS  OF  MOTION 


A.  Mean  Flow  Equations 


The  starting  point  in  the  derivation  of  the  turbulent  flow  equations  of 
motion  are  the  compressible,  laminar  Navier-Stokes  equations  (including 
continuity,  energy  and  equation  of  state). 


Conti nuity 


*  4- (,Uj)  =  ° 

J 


Momentum 


3  (pu  • )  3  *P  3t  i . 

Tt— *  THT  (puiuj)  "  Tx"  +  “ax' 

J  J  J 


Energy 


i^iii  +  *  (pHuj )  =  +  77:  (uiTij  " 


Constitutive  equations 


r  3ui  ,  3uj  z  3uk  .  i 

’j '  “  K  ^7'*^  *«J 


n  _  k  3T 
qj  '  '  k  7x 


Definition  of  total  enthalpy 


■  eint  *  P/e  *  7  ui  ui 


Sahu,  J't  "Navier-Stokes  Computational  Study  of  Axisyrmetric  Transonic 

Turbulent  Flouts  with  a  Two-Equation  Model  of  Turbulence Ph.  D . 
Dissertation j  Mechanical  and  Aerospace  Engineering  Department ,  Univer¬ 
sity  of  Delaware ,  June  1984. 
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Equations  of  state 


p  =  pRT 

(7) 

h  -  CpT. 

(8) 

Since  these  equations  cannot  be  solved  at  this  time  for  fully  turbulent 
flow  because  the  length  and  time  scales  of  the  motion  are  too  disparate,  solu¬ 
tions  are  sought  to  altered  equations  produced  by  averaging  over  time  and 
neglecting  or  modeling  various  terms.  The  derivation  used  here  is  based  on 
mass-weighted  time-averaging. 15  16  By  definition  the  mass-weighted  mean 
(denoted  by  tilde)  is 


a  =  pa/  p  (9) 


where  an  over  bar  is  a  conventional  time  average.  We  will  use  a  single  prime 
to  indicate  the  fluctuating  part  in  the  mass-weighted  variable,  i.e. : 

a(r,t)  =  a(r)  +  a'  (r,t)  (10) 


and  a  double  prime  is  used  to  denote  the  fluctuating  part  in  a  simple  time 
average: 


a(r,t)  *  a(r)  +  a"  (r,t)  (11) 


It  follows  then  that 


pa'  *  0  and  a"  =  0 


but  ~a*  *  0.  Based  on  these  definitions  one  can  derive  the  following 
equations: 


IS.  Cebeci,  T.,  and  Smith,  A.M.O.,  n Analysis  of  Turbulent  Boundary  Layers 
Academia  Press,  Inc.,  1974. 

26 .  Subesin,  M.W.,  M A  One-Equation  Model  of  Turbulence  for  Use  with  the  Cotr.- 
preseible  Navier-Stokes  Equations, n  NASA  TMX-73,  1976. 

10 
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3p  3  - 

-  -3x7  +  TxT  (?ij  ‘  puiuj) 
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(13) 


where  -"puTUT  is  the  Reynolds  stress  tensor.  These  additional  unkowns  render 

the  set  of  equations  insufficient  and  therefore  a  generalization  of  the 
Boussinesq  formula  is  introduced.  If  we  write 


•  puiuj  =  yt 


+  A  6.  . 


where  A  is  an  invariant  of  the  Reynolds  stress  tensor.  The  trace  of  the 
tensor  gives, 


A  = 


Thus,  since 


SUj 

“ij  axj  9x(  ?  9xk  ij 

we  can  write  the  Reynolds  stress  as 


2  - 


pu;u;  *  ut  dij  -  i  5  K  sij 


where  K  is  the  mass  weighted  average  kinetic  energy. 


(14) 


K  =  0.5  pu?u^  /  p  . 


The  molecular  stress  can  also  be  written  in  terms  of  d^j  as  follows. 


ij 


=  v  d 


ij 


+  ud' . 
ij 


(15) 


The  second  term  on  RHS  of  this  equation  is  the  interaction  of  the  variable 
transport  property  and  the  fluctuating  velocities  and  is  usually  neglected. 


If  this  tern  can  be  ignored  then  the  resulting  equation  for  the  combined  mole¬ 
cular  and  turbulent  stress  becomes 


?ij  -  pu;uj  ■  +  *t>  - 1 5  ^  sij  •  (*6) 

In  most  algebraic  turbulence  models  which  do  not  compute  fluctuating  quan¬ 
tities,  it  is  not  possible  to  account  for  the  last  term  involving  the  tur¬ 
bulent  kinetic  energy.  Neglecting  this  term  has  the  advantage  that  the  total 
shear  stress  is  proportional  to  the  mass  weighted  average  of  d^j  and  the  equa¬ 
tions  derived  for  compressible  laminar  flow  may  be  used  to  compute  turbulent 
flow  by  simply  replacing  the  viscosity  coefficient  u  by  the  sum  of  the  mole¬ 
cular  and  turbulent  viscosity  ii  +  u^.  This  assumption  has  been  made  in  the 

computations  reported  here  but  the  problem  is  pointed  out  and  needs  further 
i nvestigati on. 

The  turbulent  "eddy"  viscosity,  u^,  can  be  dimensionally  described  in 

terms  of  a  density  times  a  characteristic  velocity  times  a  characteristic 
length.  In  the  algebraic  turbulence  models,  the  characteristic  velocity  is 
generally  proportional  to  the  mean  velocity  gradient  so  that 


»t  •  S  i2  . 

l  dy 


(17) 


In  the  "two  equation"  turbulence  models  considered  here,  the  velocity  scale  is 
related  to  the  square  root  of  the  kinetic  energy  of  turbulence,  l!1^2.  The 
length  scale  is  equal  to  K3,/2/E,  where  E  is  the  dissipation  rate,  so  that 


-t  =  c  » K^r 


(18) 


where  c^  is  assigned  the  value  cy  =  0.09.  A  comprehensive  discussion  of  alge¬ 
braic,  one-  and  two-equation  models  of  turbulence  is  contained  in  the  survey 
by  Launder  and  Spalding.7 

The  energy  equation  is  modified  recognizing  that 


pH  =  e  +  p 


(19) 


where 


?  «v 

ujy 


!K 


■X'l 


or  the  total  energy  per  unit  volume  of  the  flow  (eint  =  internal  energy /mass  = 
cvT).  Consistent  with  the  averaging  process  used  previously,  the  turbulent 
flow  energy  equation  becomes 


-,1  (i  ♦  P)  +  -5J-  C(i  ♦  p)UjD  ■ 


♦  T^T  (-5j  -  <«  *  p)“j  *  V  dij 


+  U.  ud'  .  +  u't.  .  ), 
1  1 J  1  1J  ' 


Fourier's  law  becomes 


Q,-  =  -(k  +  k")  (T  +  r-). 

J  dXj 


The  mean  molecular  heat  flux  and  stress  work  are  taken  to  be 


.  a  ;  ,  ~  -  3 

dxT  uij  • 

J 


In  analogy  with  the  Boussinesq  representation  for  the  Reynolds  stresses  we 
write  the  following  for  the  heat  flux  and  shear  work: 


(R  +  k.  )  V—  +  jj.  t.  .  =  -5 .  -  (e  +  p)uT  +  u-t.  .  +  u't.  . 
t;  SxT  ij  J  i  ij  1  ij 


so  the  the  turbulent  heat  flux  is  determined  by 


kt  IT  =  k"  WT  •  (e  +  p)uj  +  uiTij  +  Vdi  j 

J  J 


It  is  convenient  to  eliminate  temperature  from  the  final  equation  using 


T  s  s!/yR  •  rrnr; 


which  is  averaged  to  give: 


T  -  a2/[(y-l)Cp]. 


The  thermal  conductivities  are  also  eliminated  by  introducing  the  molecular 
and  turbulent  Prandtl  numbers 


Pr  =  Cp u/k  , 


Pr^  =  CpUt/k^. 

Thus  the  final  form  of  the  energy  equation  for  turbulent  flow  becomes: 


3  ~  -n  3  1  f  5  A  wt  3a2 

TJ-  CujCe+p)]  =-5^Tf^T[1^+r^J-5J-  • 

•  t 


The  reason  for  developing  the  mean  flow  Navier-Stokes  equations  in  the 
above  form  is  to  define  and  emphasize  clearly  some  of  the  underlying 
assumptions  which  have  been  made  in  using  the  equations  in  this  form. 

1.  The  velocity  components,  u.,  are  mass-weighted  time  averages.  Since 

J 

most  solutions  do  not  provide  information^  regarding  the  density  fluctu¬ 
ations  it  is  not  possible  to  relate  u,  to  u,,  the  physical  velocity. 

J  J 

2.  A  term  involving  the  kinetic  energy  of  turbulence  has  been  neglected 
in  formulating  the  Boussinesq  approximation  to  the  Reynolds  stresses, 
i.e. : 

-  ^  p  K  S... 

Although  this  term  applies  to  the  normal  stresses  and  therefore  nay  be 
expected  to  be  relatively  unimportant,  it  should  be  possible  in  the  K-E 
models  to  test  its  significance. 

3.  y.  d..  =  yd*.  -  pu*u* 

t  ij  ij  i  j 

As  already  pointed  out  in  the  mass-weighted  average  derivation,  *  0 

*  J 

and  thus  even  if  y  =  y,  the  1THTT  tern  is  not  obviously  zero  and  y.3.. 

•  J  ^  1  j 

represents  more  than  the  conventional  form  of  the  Reynolds  stresses 

4.  The  turbulent  heat  flux,  equation  (25),  also  contains  a  number  of 
additional  terms  which  would  not  appear  in  incompressible  flow  even  if 
k"  is  negligible. 
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5.  The  definition  of  the  mean  total  energy  can  be  used  to  provide  an 
equation  for  the  pear  pressure 

p  =  (y-1  )  [  e  -  j  3ujui  -  y  puTuI  ]  . 

In  most  computational  work,  the  kinetic  energy  of  turbulence  in  this 
equation  is  ignored  in  the  evaluation  of  the  mean  pressure.  This  is 
necessary  when  using  the  algebraic  turbulence  models  but  could  be 
accounted  for  when  using  the  K-E  equations. 

6.  In  general,  the  fluctuations  in  transport  properties  k"  and  M"  are 
neglected  (see  for  example  References  17  and  18).  However,  the  only 
justification  is  mathematical  convenience. 


Turbulent  Kinetic  Enerc 


luation 


Multiplying  the  momentum  equation  by  the  instantaneous  velocity,  taking 
the  time  average,  considering  only  the  terms  where  j  =  i  and  subtracting  the 
mean  kinetic  energy,  one  obtains  an  equation  for  the  kinetic  energy  of  the 
turbulent  fluctuations 


The  four  terms  on  the  right  hand  side  represent  the  diffusion  of  turbulent 
kinetic  energy,  tendency  toward  isotropy,  dissipation  and  production  terms 


17,  White  ,  F.M. ,  "Viscous  Fluid  Flow/'  McGraw-Hill  Book  Company,  1974, 

18.  Schubauer,  C.B.,  and  Tchen,  C.M. ,  "Turbulent  Flaw,"  Turbulent  Flow  and  Heat 
Transfer,  Ed.  C.  C.  Lin,  Princeton  Series  on  High  Speed  Aerodynamics  and 
Jet  Propulsion,  Princeton  University  Press ,  1959. 


respectively.  In  this  derivation,  the  dependent  variable  is  the  kinetic 
energy  per  unit  volume  pulul/2  obtained  for  incompressible  flow.  In  the  limit 
that  p  -*■  constant  and  thus  ul  -*■  u'.'  with  t.  .  =  +  t!'  the  above  equation 

'  J  *  J  *  J 

reduces  to  the  usual  incompressible  form.  Thus  for  the  compressible  case  it 
is  only  necessary  to  rewrite  the  kinetic  energy  in  terms  of  the  mass  weighted 
average 

K  =  0.5  pujut/p  .  (33) 

The  tendency-toward- isotropy  term  is  ignored  in  forming  a  model  equation  to 
represent  equation  (32).  The  dissipation  is  also  considered  as  a  mass- 
weighted  average  so  that 

8ui  -  ~ 

T1j  Tx1  ■  p  £  <34> 

j 

and  finally  the  production  term  is  written  using  the  Boussinesq  form  as 

(35) 

(the  contribution  of  u"dl .  has  also  been  neglected).  Thus  the  model  kinetic 
energy  equation  becomes 

T5t  (p  K)  =  73T  [  [5  +  o^]  777  ]  '  3  e  +  wt  ^ik  Txj;  (36) 

where  according  to  Jones  and  Launder10,  =  1.0. 

C.  Dissipation  of  Turbulent  Energy  Model  Equation 

The  dissipation  rate  of  turbulent  kinetic  energy  e  is  written  in  analogy 
to  the  model  kinetic  energy  of  turbulence  equation  as: 

T5t(5  e)  =  ^[l5  +  "T  ]  7x7  ]  ‘  c25  =“  +  C1  |  vAk  73^  *  (3?) 

The  production  and  dissipation  terms  have  the  same  form  as  their  kinetic 

energy  equation  counterparts  but  scaled  by  */K  and  multiplied  by  constants:10 
Cj  »  1.44,  C2  *  1.92  and  oe  =  1.3. 


D.  Low  Reynolds  Number  Effects 


In  the  form  indicated  by  equations  (36)  and  (37),  the  dissipation  remains 
finite  at  the  wall  and  is  balanced  by  the  diffusion  term.  It  is  desirable  for 

computational  purposes  to  redefine  e  such  that  its  value  at  a  bounding  wall  is 

zero.  This  is  achieved  by  adding  a  term  to  the  definition  of  e  which  is 

according  to  Jones  and  Launder10 


If  e  is  replaced  by  t  in  the  dissipation  rate  equation,  the  diffusion  term 
remains  finite  at  the  wall.  The  equation  can  be  balanced  by  adding  a  tern  as 
fol 1 ows 


The  modified  K-E  equations  becomes 


Van  Gulick11  using  a  series  expansion  technique  has  shown 

that  K  =  0(y2)  as  y  *  0,  and  that  £  *  0(y)  with  the  result  that  vt  =  0(y3)  as 

suggested  by  White. 17  However,  van  Gulick' s  experience  with  the  boundary 

layer  equations  showed  that  the  K1^2  term  made  the  numerical  solution  highly 
unstable.  As  a  result,  he  adopted  a  proposal  by  Chien19  in  which  c  was 
changed  to 


1,9.  Chien,  K-Y.,  "Predictions  of  Channel  and  Boundar\  Layer  Ficus  with  a  Lou - 
Reynolds-N umber  Turbulence  Kodel,"  AIAA  Journal 3  Vol .  20,  No.  1,  January 
1982. 


t  =  e  -  2v  K/y2 


(41) 


and  the  added  term  in  the  E  equation  becomes 

-2v  i-  exp  (-  i  /)  (42) 

>2  2 

where  y+  =  u_y/v  .  This  term  becomes  small  very  rapidly  as  y  increases  and 

1  W 

has  the  correct  form  to  balance  the  diffusion  term  at  y  s  0.  Chien  also 
suggested  that  cy  in  the  equation  for  ut  include  a  damping  term  similar  to 

that  of  the  van  Driest  damping  factor 


cy  =  0.09  (1-  exp  (-0.01  y+) ) 
which  maintains  =  0(y3)  near  y  =  0. 


E.  Rotta's  Three-Dimensional  Stress  Tensor  Model 


The  Boussinesq  postulation  used  to  simplify  the  writing  of  the  Reynolds 
stress  terms  in  the  momentum  equation  implies  that  the  eddy  viscosity  produced 
by  turbulence  is  isotropic  and  that  the  differences  in  the  various  components 
of  the  Reynolds  stress  are  just  determined  by  the  components  of  the  strain 
rate  tensor,  d-j j .  In  general,  the  full  equation  for  the  Reynolds  stress  is  a 

complicated  tensor  relation  which  could  permit  Reynolds  stresses  in  three 
dimensions  which  are  not  aligned  with  the  local  mean  rate  of  strain. 

Employing  the  incompressible  equations,  Rotta20  asserts  that  the  non¬ 
isotropy  of  the  flow  is  due  to  the  pressure  strain  term  which  he  writes  as  two 
separate  components: 


p"  r  du"  dv"  -I  _  pm~f  du"  ay"] 

p  1 3y  Sx  J  "  p  L  *x“'  J 


]• 


20.  Rottc,  J.  C.  ,  "A  Family  of  Turbulence  Models  for  Three-Dimensionc l  Thin - 
Shear  Layers,"  Symposium  of  Turbulent  Shear  Flovs,  University  Park, 
Pennsylvania,  1977. 
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The  terms  involving  p 


were  assumed  to  be  proportional  to  the  mean  strain 


rates  in  streamline  coordinates  and  to  \"z.  After  transforming  from 
streamline  to  general  Cartesian  coordinates  he  found 


where 

Xxx  =  ^ss°2  +  *nn*2)/(°2  +  ^ 

Xxz  =  {xss  *  xnn}  D  +  *2)  <45) 

xzz  =  (xss^  +  xnnC2)/(Cl2  + 


and  where  x$s  and  xnn  are  constants  of  proportionality  in  the  streamline 
coordinate  system. 


The  other  factors  in  the  pressure-strain  terms  involving  p  ''  were 
presumed  to  be  proportioned  to  the  local  Reynolds  stress,  i.e. : 


pt"  r  3u"  .  3v"i 
p  Uy  ax  J 


u"v" 


/l 


pt"  r  aw"  av"i 

p  L  ay  +  "SI  J 


/L 


(46) 


where  Kp  is  a  constant  of  proportional ity.  Neglecting  all  the  other  terms  in 

the  equation  for  the  Reynolds  stresses  Rotta,  determined  the  stress  components 
as: 


-pu"v" 


(uZ+w^) 


-pw"v" 


(  tu2  *  F-w2)  *  (1-F)  u  i  ||] 

[(1-F)  0  *  (w2  ♦  F  C2)  4|] 
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(47) 


mssmmmmmjmwmaiiaw 


where  p^  is  the  isotropic  eddy  viscosity  corresponding  to  F  =  1  and  where,  in 
general , 


1  -  x. 


F  is  thus  a  parameter  to  he  selected  by  comparison  with  experimental  data. 
Rotta  suggested  a  value  of  0.5. 

F.  Effect  of  Non-Isotropy  on  the  K-Equation  Production  Term 


Van  Gunck  applied  the  Rotta  theory  to  the  boundary  layer  form  of  the  K-E 
equations  where  the  production  term  can  be  written  as 

P  =  PU'V'  +  PW'V'  .  (49) 

If  it  is  assumed  that  the  above  analysis  holds  at  least  approximately  for  the 
compressible  case,  then  P  can  be  written 


Isotropic  part 


Non- Isotropic  part 


In  the  limit  as  y  +  0,  u  ♦  0  and  w  *  Or  (circumferential  velocity  due  to 
spin) : 


p  ♦  p 


This  shows  that  near  y  =  0  where  may  be  large,  the  production  term  is 

strongly  influenced  by  the  non- isotropy.  But  experiments  show  that  w,  the 
cross  flow  velocity  decreases  very  rapidly  away  from  a  spinning  surface  and 
for  y+  >  5 

w2  -  0.1  u2  andl^l2  -  0.01  f-^12  . 


so  that  for  y+  >  5 


P  =  "t  [l|f  <F[|f]  •  <«> 

or  essentially  isotropic  for  the  boundary  layer  situation.  Because  of  the 
rapid  fall  off  of  the  effect  of  spin  on  the  cross  flow  velocity  and  the  small 
cross  flow  velocities  at  small  angles  of  attack,  the  non-isotropic  effects  can 
be  expected  to  be  very  slight  even  for  values  of  F  significantly  different 
from  unity. 


G.  Algebraic  Eddy  Viscosity  Models 


Two  forms  of  algebraic  eddy  viscosity  were  employed  for  comparison  with 
the  two-equation  K-E  model.  For  the  boundary  layer  equations,  a  form  of  the 
Prandtl-van  Driest  mixing  length  approach  was  used 


(53) 


where 


1/6  - 1  tanh[ii][1  -  exp  [-£•]] , 

X  =  0.09,  <  =  0.4,  A+  =  26. 


(54) 


The  other  model  is  the  more  sophisticated  Baldwin-Lomax21  technique,  which  was 
used  in  connection  with  the  PNS  code.  This  two  layer  model  is  based  on  the 
proportionality  between  u  and  the  local  vorticity  in  the  "innner"  region, 
such  that 


(uj.  =  p  l2 
'  t'mner 


(55) 


where 


*  »  «y  Cl  -  exp(-y+/A+)], 


(56) 


21.  Ealcbin,  B.S. ,  and  Lomax t  E. ,  ” Thin  Layer  Approximation  and  Algebraic  Mode 
for  Separated  Turbulent  Flove,"  AIAA  Paper  78-257 ,  January  1978. 
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B.  Boundary  Layer  Results 


The  results  obtained  with  the  Dwyer  and  Sanders  boundary  layer  code 
incorporating  the  K-E  turbulence  model  are  fully  described  in  van  Gulick's 
thesis11  and  in  a  paper  by  van  Gulick  and  Danberg. 12  Their  results  are  only 
summarized  here  in  regard  to  the  comparison  with  the  experimental  data  of 
Kayser  and  Sturek22  and  in  regard  to  the  effects  of  non- isotropy  and  buoyancy. 


1.  Experimental  Data 

Figure  3  shows  the  configuration  tested  in  the  Ballistic  Research 
Laboratory  wind  tunnel  by  Kayser  and  Sturek.  The  tests  were  performed  at  M  = 
3.0,  Reynolds  number  of  9.68  x  106  per  meter  and  angles  of  attack  of  0°  and 
2°.  Measurements  were  made  with  and  without  the  model 

spinning  (fl  =  20,000  RPM).  Wall  pressure  distributions  were  obtained  for  the 
non-spinning  case  and  impact  probe  boundary  layer  surveys  provided  the  major 
experimental  results.  The  impact  probe  data  were  reduced  to  provide 
longitudinal  velocity  profiles  assuming  constant  pressure  across  the  boundary 
layer  and  employing  a  Crocco  temperature-vel ocity  relationship  in  lieu  of 
measured  temperature  profiles.  A  number  of  stations  were  investigated  but 
results  will  be  limited  to  just  one  station  at  3.33  caliber  from  the  nose. 
This  station  is  on  the  cylindrical  mid-body  just  aft  of  the  ogive-cylinder 
juncti on. 


2.  Comparison  with  Experimental  Velocity  Profiles 

Figure  4  shows  a  typical  set  of  velocity  profiles  at  the  X/D  = 
3.33  station.  The  experimental  data  are  compared  with  the  computed  three- 
dimensional  boundary  layer  code  predictions  using  the  algebraic  and  K-E 
turbulence  models.  The  comparison  is  based  on  the  physical  variables  and 
therefore  demonstrates  the  ability  of  the  computational  system  (inviscid  and 
viscous  flow)  to  predict  the  measured  quantities.  Displacement  effects  are, 
however,  not  included.  The  predicted  velocities  are  generally  less  than  those 
measured.  The  K-E  prediction  is  in  somewhat  better  agreement  with  the  data 
than  that  of  the  algebraic  model.  The  computed  boundary  layer  thickness  is 
larger  than  the  experimental  data  which  is  in  part  caused  by  the  incorrectly 
calculated  pressure  history  obtained  from  the  theoretical  prediction. 
At  <|>  =  180°,  the  lee  side,  the  agreement  is  better  at  all  stations  investi¬ 
gated  and  the  lee  side  experimental  data  are  expected  to  be  more  reliable 
because  of  the  thickness  of  the  profile. 

A  comparison  is  also  shown  in  Figure  5  between  the  experimental  profile 
data  and  prediction  based  on  non-dimensional ized  variables.  The  local 
velocity  is  normalized  by  the  boundary  layer  edge  velocity  and  the  normal 


22. 


Kauser ,  L.D. ,  and  Sturek ,  W.B. ,  "Experimental  Measurements  in  tte  Turbulent 
Boundary  Layer  of  a  Yawed,  Spinning  Ogive-Cylinder  Body  of  Revolution  at 
Mach  3.0.  Part  II:  Data  Tabulation US  Army  Ballistic  Research 
Laboratory ,  Memorandum  Report  ARBRL-MR-02813 ,  March  1978. 
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distance  by  the  boundary  layer  edge  velocity  and  the  normal  distance  by  the 
boundary  layer  thickness.  This  presentation  concentrates  on  the  ability  of 
the  prediction  to  reproduce  the  form  of  the  turbulent  profile  independent  of 
any  discrepancy  in  the  outer  boundary  condition.  The  results  show  that  the 
K-E  method  is  generally  superior  to  the  algebraic  solution  although  the 
algebraic  solution  is  at  most  7  percent  lower  than  the  K-E  solution. 

Figures  6  and  7  show  typical  K  and  E  profiles  close  to  the  wall  obtained 
from  these  calculations.  These  results  are  non-dimensional ized  in  terms  of 
1 aw-of-the-wall  variables 


(K+  =  K/u*,  E+  =  v^/uj,  y+  =  uTy/V  where  u t  =  Vp>1/2>* 


3.  Effect  of  Non- Isotropy 

Rotta's  non-isotropy  formulation  of  the  turbulent  stresses  was 
incorporated  into  the  boundar7  layer  mean  flow  equations.  The  K-E  or 
algebraic  turbulence  model  provided  the  which  is  multiplied  by  an  effective 

velocity  gradient  composed  of  an  isotropic  and  non- isotropic  part 


3u 


-  (1-F) 


t  cly  J  *  Ty 

effective 

Isotropic 


~r~  3u  ~  3w7 

3U  T  "Sy  w'5yJ 

Ty - p? - 

y  (u  +  w  ) 

Non- Isotropic 


(65) 


Figure  8  shows  a  plot  of  the  isotropic  and  non-isotropic  terms  for  a  typical 
solution.  The  non-isotropic  part  is  seen  to  fall  off  dramatically.  The  net 
effect  on  the  longitudinal  flow  is  found  to  be  small  for  all  the  conditions 
considered.  In  varying  the  constant  F  between  1  and  0.25,  the  u- velocity 
profiles  and  longitudinal  Cfx  do  not  change  discernably  (see  Figure  9). 

However,  the  circumferential  flow  is  significantly  affected  as  shown  in  Figure 
10  where  the  circumferential  skin  friction  coefficient  is  plotted  against 
angular  position  for  the  spinning  body.  The  general  effect  of  decreasing  F  is 
to  increase  nearly  linearly  and  uniformly  the  circumferential  friction  (in 
this  case  to  decrease  the  calculated  roll  damping).  The  longitudinal  skin 
friction  decreases  slightly  because  of  this  formulation  of  the  shear 
stresses.  This  implies  an  exchange  of  momentum  from  longitudinal  to  the 
circumferential  direction  with  decreasing  F  consistent  with  the  assumptions 
underlying  Rotta's  stress  tensor. 


4.  Effects  of  Streanline  Curvature 


Bradshaw23  has  proposed  an  effect  of  streanline  curvature  on 
turbulence  based  on  an  analogy  with  buoyancy  effects.  Bradshaw's  proposal 
called  for  a  nodification  of  the  nixing  length  by  a  linear  factor  involving 
the  turbulent  Richardson  number. 

The  fundamental  interpretation  of  Richardson  number  is  that  of  a  ratio  of 
buoyancy  to  inertia  forces.24  For  example.  Reference  6  cites  a  formulation  by 
Prandtl,  as: 

Ri  =  "9  [p  ty  ]  /  [ly]  *  (66) 


Bradshaw  has  shown  that  this  can  be  interpreted  as  the  ratio  of  two  time 
scales  squared.  He  then  defined  a  time  scale  for  the  mean  flow  for  curved 
streamlines  as 


(mean  flow  time  scale)^  -i 
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where  V  /r  is  the  radial  acceleration  and  the  term  in  the  square 

the  inverse  of  a  length  scale  based  on  the  rate  of  change 
momentum.  He  compared  this  to  a  time  scale  for  turbulence 
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Thus  the  Richardson  number  for  turbulent  flow  along  curved  streamlines 
defi ned 


is 
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Bradshaw  also  noted  that  compressibility  would  effect  the  turbulent  Richardson 

y-1  2 

number  and  he  recommended  it  increase  in  proportion  to  (1  +  M  )  . 


23.  Bradehau',  F.,  "The  Analogy  Between  Streamline  Curvature  ana  Buoyancy  in 
Turbulent  Shear  Flov,"  J.  Fluid  Mech. ,  Vol.  36,  Part  I,  1969. 

24.  Arpaci ,  V. S. ,  ana  Larsen,  P.S.,  "Convective  heat  Transfer, •’  Prentice- Hell, 
1984. 
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Launder  et  al25  have  attempted  to  extend  these  ideas  to  the  K-E  theory  by 
redefining  the  turbulent  time  scale  as  K/E  so  that 

f  .  2  ^  q  )  ~2l  v_i  2 

Rit  =  [  <KV  '  (r  E)  J  C1  +  V  M  }  •  (68) 

It  was  suggested  that  the  most  appropriate  way  for  the  Richardson  number  to 
influence  the  Reynolds  stresses  was  through  modification  of  one  or  more  of  the 

empirical  constants  in  the  E  equation.  Specifically  it  was  recommended  that 
C2  be  increased  by  the  factor 


1  -  Cc  Rit  (69) 

where  a  value  of  Cc  =  0.2  was  cited  as  optimum  based  on  results  for  low  speed 
rotating  body  flow  fields. 

Calculations  were  made  for  the  wind  tunnel  configuration  previously 
described  at  2°  angle  of  attack  and  with  spin.  The  Richardson  number  was 
estimated  locally  in  terms  of  the  circumferential  velocity,  w,  because  it  was 
anticipated  that  the  major  effect  would  be  associated  with  the  circumferential 
flow  near  the  model  surface.  The  Ri^  employed  was 

Rit  *  2f?  7  [If]  (1  •  (70) 

This  modification  of  the  dissipation  in  the  E  equation  had  very  little 
effect.  The  longitudinal  skin  friction  was  slightly  increased  as  shown  in 
Figure  11.  The  magnitude  of  the  streamline  curvature  effect  is  negligible  for 
the  present  situation.  This  result  might  be  different  at  higher  angles  of 
attack. 


IV.  APPLICATION  TO  THE  PARABOLIZED  NAV I ER- STOKES  CODE 


A.  Mean-Flow  Equations 

The  parabolized  Navier-Stokes  method  (PNS)  employed  here  was  developed  by 
Schiff  and  Steger2  as  a  technique  for  calculating  steady,  supersonic,  high 
Reynolds  number  flow  about  three-dimensional  configurations  with  moderate 
axial  geometry  variation.  The  PNS  equations  are  obtained  from  equations  (12), 
(13)  and  (30)  which  can  be  written  in  the  following  compact  form: 


25.  Launder j  B.E. }  Priddin ,  C.  H. ,  and  Shama,  B.F.  t  '"The  Calculation  of  Turbulent 
Boundary  Layers  on  Spinning  and  Curved  Surfaces J.  Fluids  Engineering , 
pp.  231-239 ,  March  1977. 
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The  physical  x,  y,  z  Cartesian  coordinate  system  is  transformed  into  new 
body  oriented  computational  coordinates  S,  n,  c  as  illustrated  in  Figure  12. 
This  transformation  simplifies  the  application  of  boundary  conditions  and 
makes  it  possible  to  neglect  certain  viscous  terms  because  of  the  high 
Reynolds  number  -  thin  shear  layer  approximation. 2  26  The  transformation 
relations  are: 


5  =  £(x) 
n  =  n(x,y,z) 
C  =  S(x,y,z) 
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The  resulting  PNS  equations  become: 
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26.  Steger,  J.L.,  ", Implicit  Finite-Difference  Simulation  of  Flow  About 
Arbitrary  Two-Dimensional  Geometries ,"  AIAA  Journal,  Vol.  16,  No.  7  , 
July  1978.  ~~~  "  ~—™° 
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These  equations  have  the  strong  conservation  form  while  retaining  the 
Cartesian  velocity  components  as  dependent  variables.  The  contravariant 
velocity  components  U,  V,  W  are  defined  as: 


U  =  C  u 
x 


V  =  o  u  +  n  v  +  n  w 
x  y  z 


W  =  ;  u  +  c  V  +  s  w  . 
x  y  z 


(73) 


The  thin  viscous  layer  approximation  retains  only  those  terms  involving 
derivatives  with  respect  to  c,  the  near  normal  to  the  wall  coordinate.  Thus 
the  RHS  terms  can  be  written  as: 
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jgration  scheme  used  i 
factorized  algorithm  in  delta  form  developed  by  Beam  and  Warming.3  27  The 


The  numerical  integration  scheme  used  is  the  implicit  approximate 

i  nn  .  3 


27 ,  Bean,  R.M.,  and  Warming,  R.F. ,  "An  Implicit  Factored  Scheme  for  tr.e 
Compressible  Navier-Stokes  Equations,"  AIAA  J.,  Vol,  16,  No.  4,  April 
1978. 


difference  equations  are  linearized  in  their  vector  form  and  the  solution 
requires  block  tridiagonal  inversion  at  each  marching  step.  See  reference  28 
for  additional  details  and  applications. 


B.  Turbulence  Models 


The  Baldwin-Lomax  algebraic  model21  was  incorporated  into  the  method  as 
it  appears  in  the  work  of  Schiff  and  Sturek.29  The  K-E  equations  (38)  and 
(39)  were  incorporated  by  Kim13  in  the  steady,  supersonic  P NS  code  paralleling 

the  work  of  Sahu14  in  his  application  of  the  K-E  equations  to  the  time- 
dependent,  transonic  Navier-Stokes  problem.  A  new  formulation  of  the 
turbulent  viscosity  subroutine  from  that  employed  by  van  Gulick  was  required 
because  of  the  body  oriented  coordinate  transformations  in  the  PNS  solution. 

The  K-?  equations  can  be  written  in  a  transformed  form  similar  to  the  mean 
flow  equations,  i.e. : 
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28.  Sturek t  W.B.,  and  Schiff t  L.B. ,  "Computation  of  the  Magnus  Effect  for 
Slender  Bodies  in  Supersonic  Flex),"  AIAA  Paper  80-1S86CP ,  August  1980. 

29.  Schiff t  L.E.,  and  Sturek t  W.B.,  "Numerical  Simulation  of  Steady  Supersonic 
Flou  Over  an  0  give-Cylinder-Boattail  Body  ”  AIAA  Paper  80-0066 ,  January 
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Note:  these  equations  have  been  non-dimensional  ized  using  a.,  p„,  um  (f ree- 

strean  properties)  and  the  body  reference  dianeter  D.  Thus  the  Reynolds  num¬ 
ber  parameter  is: 

Re  =  P.aj)/M.  .  (75) 

In  addition,  ?  has  been  eliminated  from  the  K  equation  and  K  from  the 
E-equation  using  the  following  definition  of  the  turbulent  viscosity 


ut  =  Re  c  p  K^/E  . 
t  u 


Thus  the  equations  are  decoupled  but  require  an  estimate  of 
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1.  Finite  Difference  Approximation 


The  governing  K-E  equations  have  a  form  similar  to  that  of  the  mean 
flow  PNS  equations.  In  fact  they  could  be  solved  by  increasing  the  order  of 
the  block  tri diagonal  inversion  scheme  of  the  PNS  code  but  this  would  have 
reduced  the  flexibility  of  the  code  and  increased  the  running  time.  The  sane 
kind  o*  linearization  and  implicit  factorized  delta  form  is  employed  for 

the  K-E  equation  with  the  exception  that  in  the  decoupled  form  the  inversion 
requires  only  an  efficient  tridiagonal  algorithm.  The  boundary  conditions 
employed  at  the  body  surface  are  the  same  as  in  the  boundary  layer  case, 

i.e. ,  K  =  E  =  0.  At  the  outer  boundary  a  zero  gradient  condition  was  applied 

iE  -  if  -  o 

3;  3;  u  * 

Thus  questions  regarding  the  determination  of  the  edge  of  the  boundary  layer 

and  evaluation  of  the  inviscid  K-E  equations  along  an  inviscid  streamline  were 
avoided.  The  above  conditions  were  applied  at  an  arbitrary  point  well  away 
from  the  viscous  layer  near  the  body. 

The  initial  conditions  employed  to  start  the  calculation  are  based  on 
equating  the  production  and  dissipation  terms.  This  is  obviously  incorrect  at 

the  wall  and  at  the  outer  edge  of  the  boundary  layer.  However,  the  error  is 

expected  to  decay  as  the  flow  develops.  The  resulting  equations  for  initial 

values  of  K  and  E  are 
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The  turbulent  viscosity  and  initial  velocity  distribution  were  obtained  using 

the  algebraic  model  ahead  of  the  station  where  the  K-E  equations  were  first 
introduced.  A  more  complete  description  of  the  finite  difference  procedure 
can  be  found  in  Kim's  Thesis. 13 


C.  Preliminary  Results  Using  K-E  Model  in  the  PNS  Code 

The  experimental  data  from  the  ogive-cylinder-boattail  projectile 
configuration  were  used  to  evaluate  the  PNS  Code.  The  calculations,  however, 
were  limited  to  zero  angle  of  attack  and  zero  spin.  The  code  was  run  on  the 
University  of  Delaware  computer  system  and  required  45  minutes  of  CPU  time  to 
complete  just  160  steps.  Because  of  the  long  running  time  on  this  system  only 
a  few  runs  were  made  in  an  attempt  to  verify  the  program. 


Figure  13  shows  one  of  the  velocity  profiles  obtained  at  the  X/D  =  3.33 
station.  The  difference  between  the  K-E  and  algebraic  computation  is  not 
significant  at  this  point.  The  agreement  with  the  experimental  data  is  not 
significantly  different  from  that  of  the  boundary  layer  solution  despite  the 
fact  that  the  PNS  calculation  included  the  inviscid  as  well  as  the  viscous 

flow  field.  The  close  agreement  between  the  two  turbulent  viscosity  models  is 
shown  in  Figure  14  where  both  methods  give  essentially  identical  results  out 

to  a  y+  of  100  near  the  maximum  in  Beyond  that  point,  the  K-E  turbulent 

viscosity  falls  off  more  rapidly  than  the  algebraic  theory  but  the  difference 

has  little  effect  on  the  computed  velocity  profile.  The  K  and  t  profiles  on 

which  the  values  are  based  are  shown  in  Figures  15  and  16.  These  show 

unexpected  sharp  peaks  and  lack  of  profile  structure.  It  is  believed  that 
further  computational  experimentation  is  needed  to  provide  additional 

substantiation  and  extend  these  results. 


V.  CONCLUDING  REMARKS 

A  two- equation  model  of  turbulence  based  on  equations  for  the  turbulent 
kinetic  energy  and  energy  dissipation  has  been  formulated  and  applied  to  the 
solution  of  supersonic,  three-dimensional  flow  and  the  results  compared  to 
experimental  data.  Two  solution  techniques  were  studied;  the  boundary  layer 
theory  approach  and  the  parabolized  Navier-Stokes  method. 

The  boundary  layer  calculations  were  made  using  the  K-E  turbulence  model 
and  compared  to  results  employing  an  algebraic  turbulent  viscosity  model. 
Both  turbulence  models  gave  good  agreement  with  velocity  profiles  obtained  on 
a  spinning  ogi ve-cyl  i nder-boattai 1  configuration.  The  K-E  model  gave  slightly 
better  results  regarding  the  shape  of  the  profile  but  both  numerical 
predictions  over-estimated  the  boundary  layer  thickness. 

The  parabolized  Navier-Stokes  technique  has  the  advantage  of  computing 
the  entire  flow  field  at  one  time  including  the  viscous  and  inviscid 
regions.  In  order  to  make  the  K-E  turbulence  model  compatible  with  the  PNS 
method,  the  turbulence  equation  had  to  be  reformulated  in  body  oriented 
curvilinear  coordinates.  Preliminary  results  with  this  code  also  showed  good 
agreement  when  compared  to  an  algebraic  turbulence  model  and  with  experimental 
velocity  profile  data. 

The  fundamental  assumptions  inherent  in  the  compressible,  turbulent  flow 
equations  were  reviewed.  The  formulation  is  based  on  mass-weighted  time 
averaged  variables,  although  no  distinction  was  made  between  mass-weighted 
time  average  and  the  measured  average  velocity  of  the  experimental  data.  In 
developing  the  final  form  of  the  conservation  equations  it  is  necessary  to 
neglect  fluctuations  in  molecular  thermal  conductivity  and  viscosity  and  other 
terms  in  the  momentum  and  energy  equations  that  arise  because  of  the  kind  of 
averaging  employed.  The  effects  of  these  assumptions  are  hidden  within  the 
subsequent  modeling  approximations  made  to  achieve  closure.  Progress  in 
developing  a  fully  consistent  system  of  equations  will  require  consideration 
of  the  importance  of  the  neglected  terms. 


An  assumption  that  has  also  been  made  in  the  compressible  flow  equations, 
is  that  the  total  energy  consists  of  only  the  thermodynamic  internal  energy 
and  kinetic  energy  of  the  mean  flow.  The  kinetic  energy  of  the  turbulent 
notion  nay  be  expected  to  be  important  close  to  any  surface  and  thus  have  an 
important  effect  on  the  mean  pressure  as  well  as  temperature  in  that  region. 

The  effect  of  incorporating  Rotta's  non- isotropic  theory  for  the  Reynolds 
stresses  was  investigated  using  the  boundary  layer  formulation.  Even  for 
widely  varying  values  of  the  independent  parameter  of  that  theory,  relatively 
little  effect  on  the  final  solution  was  observed.  The  cross  flow 
variable  c^  and  cross  flow  velocity  were  most  significantly  effected 

indicating  an  exchange  of  momentum  from  the  longitudinal  to  the 

circumferential  flow.  The  test  case  was  for  small  angles  of  attack  with  small 
non-isotropic  Reynolds  stresses  and  the  correspond- ng  effect  on  the 
dissipation  term  in  the  K-E  equation  may  be  important. 

Bradshaw's  streamline  curvature  theory  was  also  studied  in  the  boundary 
layer.  The  Richardson  number  effect  on  the  dissipation  equation  had 

essentially  no  effect  on  the  overall  solution  for  the  spinning  body  in  the 
supersonic  flow  investigated. 


Figure  1.  Boundary  Layer  Coordinate  System 
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Figure  6.  K  versus  y  of  the  Boundary  Layer  Solution 
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LIST  OF  SYMBOLS 


Speed  of  Sound 
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